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A simple theory of the fluid state of a charged colloidal suspension is proposed. The full free 
energy of a polyelectrolyte solution is calculated. It is found that the counterions condense onto the 
polyions forming clusters composed of one polyion and n counterions. The distribution of cluster 
sizes is determined explicitly. In agreement with the current experimental and Monte Carlo results, 
no liquid-gas phase separation was encountered. 
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The thermodynamic properties of systems in which 
the predominant interactions are due to the long-ranged 
Coulomb potential still remain largely not understood in 
spite of the tremendous effort that has been invested over 
the span of this century. Nevertheless, it would be un- 
fair to say that no great progress has been done. The 
pioneering work of Debye and Hiickel |^ has lead to our 
understanding of dilute electrolyte solutions. The subse- 
quent improvements by Bjerrum extended the validity of 
the limiting laws to larger densities M. These develop- 
ments were followed by the introduction of powerful inte- 
gral equations and by the computational methods such as 
Monte Carlo (MC) simulations Q . Surprisingly, the the- 
oretically obtained coexistence curve B , that is in closest 
agreement with MC, is based on the fundamental ideas 
advanced by Debye, Hiickel and Bjerrum more than 70 
years ago |l|,01. The simplicity and the transparency of 
the ideas forming the basis of the Debye-Hiickel-Bjerrum 
(DHBj) theory makes it easy to apply to other coulombic 
systems ||^. 

The charged colloidal suspensions present a severe 
challenge to any statistical mechanics theory. The asym- 
metry between the charge on a polyion and a counterion, 
which can be as high as 10000:1, makes the usual in- 
tegral equations of the liquid-state theory impossible to 
solve. For low charge asymmetry, less than 20:1, it was 
found that there exists a region in the temperature den- 
sity plane where the Hyperneted Chain equation (HNC) 
ceases to have solutions [|| . This could be interpreted as 
a region of instability, in which the sample phase sepa- 
rates into the coexisting liquid and gas. It is, however, 
still unknown to what extent the break down in HNC 
equation can be attributed to the underlying phase sep- 
aration, since the region of instability of HNC does not 
coincide with the true spinodal line W. Furthermore, the 
extensive experimental and simulation search for this gas- 
liquid transition for polyelectrolytes has, so far, proven 
to be futile §. 

At high volume fractions the strongly charged polyions 
form a lattice [hcc or fee). This "solid state" provides 
us with a major simplification in that each polyion can 
be studied individually inclosed in its own Wigner-Seitz 
cell and surrounded by its own counterions P| . Unfortu- 



nately j_once the lattice melts the cell picture is no longer 
valid [10|. As is usual, the liquid state is significantly 
more complex than the solid state. 

In this letter we shall attempt to construct a theory for 
the fiuid state of a polyelectrolyte solution. We shall work 
in the context of the Primitive Model of Polyelectrolyte 
(PMP) p]. Our system will consist of Np polyions in- 
side a volume V . The polyions will be idealized as hard 
spheres of radius a, each carrying Z ionized groups of 
charge —q uniformly spaced along the surface. A total 
of ZNp counterions will be present to preserve the over- 
all charge neutrality of the system. For simplicity, we 
shall take the counterions to be point-like and to carry 
charge -\-q. The solvent will be represented as a uniform 
medium of a dielectric constant D. As was pointed out 
by Onsager |ll|, the full nonlinear Poisson-Boltzmann 
(PB) equation is electrostatically inconsistent for the 
asymmetric systems, linearization then, besides simpli- 
fying the calculations, is an important step in restoring 
the self-consistency of the theory. The fundamental as- 
sumption behind the DHBj theory is that the nonlin- 
earities omitted in the process of linearization of the PB 
equation can be reintroduced into the theory through 
the allowance for ion association. In general, we expect 
that the fiuid state of the asymmetric electrolyte will 
be composed of free unassociated polyions of density po, 
free counterions of density pf, and of clusters consist- 
ing of one polyion and Q < n < Z associated counte- 
rions. The density of clusters with n counterions is /?„. 
In the discussion that will follow we shall suppose that 
the condensed counterion neutralizes one of the polyion 
charges, in such a way that the effective surface charge 
of a n-cluster is cr„ = ~q{Z —- n)/(47ra^). It is evident 

that pp = Y.n=oP^^ ^^d ^Pp ^ Pf + Z]„=o"''°»' where 
Pp — Np/V, is the total density of polyions (associated 
or not). All the thermodynamic properties of the poly- 
electrolyte solution can be determined once the free en- 
ergy is known. In particular the osmotic pressure is a 
Legendre transform of the Helmholtz free-energy density, 
/ = -F/V, p = f{T, {ps})+J2s P-sPs, where the chemical 
potential of a specie s G {/, "} is Ps ~ —df/dps- The 
free energy can be expressed as a sum of electrostatic 
and entropic contributions. The electrostatic free en- 
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ergy is due to the polyion-counterion, the polyion-polyion, 
and the counterion-counterion interactions. The polyion- 
counterion contribution can be obtained in the frame- 
work of the usual DH theory. 

Let us concentrate our attention on one cluster of size n 
fixed at r = 0. Due to the excluded volume, no countcri- 
ons will be found inside r < a. Therefore the electrostatic 
potential '^„ hi this region satisfies the Laplace equation 
V^^n = 0. Outside r > a — e, the mean charge distribu- 
tion will be specified by the cluster-counterion correlation 
function. Within the DH theory this is approximated by 
a Boltzmann factor leading to the charge density pq (r) — 

-j:Lo9iZ-n)pn + qPfe-^'''^"'^'-^+anS{\f\-a). Notice 
that only free unassociated counterions get polarized; the 
free unassociated polyions and clusters are too massive to 
be affected by the electrostatic fluctuations and only con- 
tribute to the neutralizing background. Substituting this 
expression into the Poisson equation, V^^„ = —Airpg/D, 
one obtains the nonlinear Poisson-Boltzmann equation. 
After the linearization of the exponential factor, we are 
lead to the Helmholtz equation, V^Vl/n = K^'^n, where 



iATTp*f/T* and the reduced temperature and den- 
sity are respectively T* = aDkBT/q^ and p* = pa^ . 
In principle the linearization is valid only in the limit 
Pq^n "^ 1, however, when the formation of clusters is 
properly taken into account the validity of the theory 
extends far into the nonlinear regime B . 

Both the Laplace and the Helmholtz equations can 
now be solved, supplemented by the boundary condi- 
tion of continuity of the electrostatic potential, and dis- 
continuity in the normal component of the electric field 
related to the presence of surface charge at r = a 
fl^ . Under these conditions, it is easily found that the 
electrostatic potential of a n-clustcr is '^^^{qn,qZ) = 



-{Z - n)q/{Da{l + na)) for r < a, and *°"'(<7n,gZ) = 
-{Z - n)qe'^'-''-'''> l{pr{\ + na)) for r > a. The electro- 
static energy of a n-cluster is 

/•oc 

(3Un{qn,qZ)^2TT(3 Pg{r)'i>n{qn,qZydr , (1) 



where P — l/{kBT). The electrostatic free-energy den- 
sity for the polyion-counterion interaction is obtained 
through the Debye charging process, where all the parti- 
cles are charged from to their final charge ||l|,Q , 
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The polyion-polyion contribution to the free energy 
can be calculated in the spirit of the usual Van-der- 
Waals (VdW) theory [Q. In particular, the tremen- 
dous difference in mass between the polyion and the 
counterion, nic/rnp <^ 1, leads to the effective separa- 
tion in characteristic time scales Tc <^ Tp. Thus, any 
change in the configuration of one polyion will be ac- 
companied by an almost instantaneous rearrangement of 
the counterion cloud. Under these conditions, the de- 
grees of freedom associated with the motion of coun- 
terions can be effectively integrated out, resulting in a 
s/iort-ranged effective pair potential of a DLVO form |14[] , 
V^%^ = q^{Z - n){Z - m)[6l(Ka)]2e-''7(L>r), where the 
enhancement factor dina) = e'"'/(f + na) is the result of 
the absence of screening inside the volume occupied by 
two polyions. Within the VdW theory 
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Now, the electrostatic free energy due to the inter- 
actions between free ions can be estimated as that of 
a One Component Plasma (OCP) ||T^. In particu- 
lar, although within our treatment the counterions are 
taken to be point particles, the electrostatic repulsion 
between the two equally charged counterions will keep 
them from approaching each other closer than a distance 
d. This can be calculated using the OCP theory to be 
d=[{l + 3Ka/T*y/^ - 1]/k [16(a)]. As one might have 
expected, in the limit of small densities this reduces to 
the Bjerrum length, d « a/T* = Xb- The electrostatic 
free energy is found through the Debye charging process 
and an analytic expression valid over a wide range of 
coupling strengths is presented in Ref. [16(b)]. 



Finally, the entropic (mixing) contribution to the free 
energy is expressed using the Flory theory as PJen = 
'^g[ps — Ps lii(0s/Cs)]j where 4>s is the volume fraction oc- 
cupied by each specie s; 0„ — 47r/9„a'^/3, 4>f = Anpfd'^/S, 
and Cs is the internal partition function oi s, (^q — Q — 1, 
and C« = [Z\/{{Z-ny.n\)]e-'^^". Here E„ is the electro- 
static energy of n counterions condensed onto the surface 
of a polyion and can be obtained through the charging 
process. 
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The minimization of the total free energy / = /en + 
fpc + fpp + fee leads to the law of mass action 



/Lto + TT'fj.f = jjLn- This is a set of Z coupled nonlinear 
algebraic equations. We were able to solve these itera- 
tively, starting with a uniform distribution of clusters. A 
sample of the distribution obtained is presented in the in- 
set of Fig. 1. The average cluster size is (n) — ^ npn/ Pp, 
and the average cluster charge is {Zcs) = Z — {n). We 
note that the width of the distribution remains quite nar- 
row and is not very sensitive to the variations in density 
or temperature. This suggests that the polydispersivity 
in cluster sizes is not very important and can be replaced 
by one characteristic cluster size. In this case the theory 
becomes extremely simple, since the free energy will be 
an explicit function of the number of bound counterions 
• • riB, since p„ — PpSn.uB- The thermodynamically sta- 
^H ble cluster, n, will be the one for which the free energy, 
^^ f{Z, Ub), attains a minimum, f{Z, n) = miun^^ f{Z, ns)- 
Thus, we must solve only one algebraic equation, instead 
of Z coupled ones. Indeed, as expected, the one-cluster 
approximation is in excellent agreement with the full the- 
ory. In Fig. 1 we present the average cluster charge as 
a function of the polyion bare charge, while in Fig. 2 
the osmotic pressure inside the polyelectrolyte solution 
is computed, which is a monotonically increasing func- 
tion of the density of polyions. 
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FIG. 1. Average cluster charge, {Z^s) as a function of the 
bare charge Z for a temperature T* = 100 and density of 
polyions Pp = 3/(327r). The solid line is calculated using 
the full theory (distribution of cluster sizes), while the circles 
were obtained within one-size cluster approximation. Inset 
depicts the distribution of cluster sizes for Z — 1000. The 
dashed line is the bare charge. 

An interesting question that arises is, what happened 
to the phase transition, so successfully predicted by the 
DHBj theory in the case of symmetric electrolyte Q]? 
The answer to this question is far from clear. However, 
our derivation does shed some light on the mechanism of 
the disappearance of the phase transition. The funda- 
mental postulate that only the counterions are polarized 
by the electric field fluctuations, can be seem to lie be- 
hind the disappearance of the transition. Indeed if this 



postulate is waved, so that both polyions and counteri- 
ons can be polarized [[l7| , it is a simple matter to show 
that the pure DH theory predicts that the suspension 
will undergo a phase separation when the temperature 
is reduced below T* — Z/16. The fundamental question 
of what is the maximum charge asymmetry above which 
the polyelectrolyte solution will remain stable still needs 
to be answered. 
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FIG. 2. Dependence of the renormalized dimensionless 
pressure Ppa''/Z on the density of macroions p* for several 
values of Z{100, 1000, 10000). The solid lines were calculated 
using the full theory (distribution) , while the circles were ob- 
tained using one-size cluster approximation. 
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